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Abstract 

The dynamics of planar crack fronts in hetergeneous media near the critical 
load for onset of crack motion are investigated both analytically and by nu- 
merical simulations. Elasticity of the solid leads to long range stress transfer 
along the crack front which is non-monotonic in time due to the elastic waves 
in the medium. In the quasistatic limit with instantaneous stress transfer, the 
crack front exhibits dynamic critical phenomenon, with a second order like 
transition from a pinned to a moving phase as the applied load is increased 
through a critical value. At criticality, the crack-front is self-affine, with a 
roughness exponent C, = 0.34 it 0.02. The dynamic exponent z is found to 
be equal to 0.74 it 0.03 and the correlation length exponent u = 1.52 it 0.02. 
These results are in good agreement with those obtained from an epsilon ex- 
pansion. Sound-travel time delays in the stress transfer do not change the 
static exponents but the dynamic exponent z becomes exactly one. Real elas- 
tic waves, however, lead to overshoots in the stresses above their eventual 
static value when one part of the crack front moves forward. Simplified mod- 
els of these stress overshoots are used to show that overshoots are relevant 
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at the depinning transition leading to a decrease in the critical load and an 
apparent jump in the velocity of the crack front directly to a non-zero value. 
In finite systems, the velocity also shows hysteretic behaviour as a function 
of the loading. These results suggest a first order like transition. Possible 
implications for real tensile cracks are discussed. 
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I. INTRODUCTION 



The dynamics of cracks in heterogeneous media is a very rich field involving much physics 
that is yet to be understood. Even in situations in which the path of a crack is predetermined 
- for example by a pre-weakend fault- its dynamics can still be complicated. The simplest 
situation is a crack confined to a plane. For small loads across such a planar crack, the crack 
front will be at rest. As the load is gradually increased, the crack front may undergo some 
transient motion but then again be arrested. If the load is increased above a critical load, 
however, the crack front will begin to propagate through the sample. The behaviour near 
to the onset of propagation of planar cracks - in particular tensile cracks - is the subject of 
this paper. 

In recent years there has been considerable theoretical progress towards understanding 
the dynamics of elastic manifolds moving through random media, such as charge density 
waves [|l|, fiuid-surface contact lines and interfaces between two phases. All of these 
exhibit a type of non-equilibrium critical phenomenon near to the onset of motion. However 
there are various features which make the system of a planar crackfront moving through a 
heterogeneous medium different from these other systems. 

For cracks (as well as for contact lines) the bulk degrees of freedom lead to effective 
long range interactions between the points on the front Thus, when a point on the 

crack front moves ahead, the stress at all other points on the front increases due to the 
elastic interactions tending to pull them forward. In addition, elastic waves are emitted as 
the crack front moves non-uniformly. When one point moves ahead, these waves result in 
stresses elsewhere on the front which, for a while, are greater than those due to just the 
static elastic deformations which will obtain long after the waves have passed. Both these 
stress overshoots and the long range interactions have earlier been shown to play a crucial 
role in the dynamics of the crack front when it is moving with a non-zero mean velocity. 
i 

In the absence of these stress overshoots — as obtains if the stress transfer is quasistatic 



3 



— many aspects of the dynamics of a planar crack front near the onset of motion can be 
understood by analogy with interfaces, in particular via a functional renormalization group 
analysis, which for cracks, as for contact lines, entails an expansion about two dimensions 
0. The phenomenology is built on the existence of two "phases" which are separated 
by a unique critical load. When the applied load, G°°, is small, there is no steady-state 
motion and the crack front is pinned by the random toughness in one of many locally stable 
configurations — we will ignore here and henceforth the effects of thermal creep. As the 
load is increased adiabatically, there are a series of local instabilities of the crack front 
which lead to "avalanches" that can become large as is increased further. Eventually at 
the critical load, G^**, the crack front de-pins and begins to move — albeit very jerkily — 
with a non-zero, mean steady state velocity, v. In an infinite system, this transition from 
the stationary to the moving phase exhibits non-equilibrium dynamic critical phenomena 
somewhat analogous to those near conventional second-order transitions. One macroscopic 
manifestation of this is the behaviour of the mean velocity of the crack front at a load just 
above the critical load: 

t; ~ - Gf )^. (1) 

A natural question that arises is the role of the stress overshoots left out of the quasistatic 
analysis. In particular, what are their effects on the crack dynamics and how do these affect 
the depinning transition? The temporal shape of the stress overshoots seen by a point on 
the crack front depends on various microscopic details, such as the microscopic response 
time of the crack front, acoustic damping processes etc. How the dynamics of the crack 
front depends on the nature of the stress overshoots and if there is any limit in which the 
dynamics of the medium can be neglected are not understood; these are questions that must 
be addressed. In particular, in the presence of stress overshoots, is there a regime in which 
a second order like transition from the pinned to the moving state persists or does the crack 
front always jump directly to a finite velocity? If the stress overshoots are "relevant" at the 
depinning transition do they make the velocity versus loading curves hysteretic? In either 

4 



case, is the 'moving phase' just above the threshold a non-trivial statistically stationary state 
or is it characterized by noisy linear dynamics? Thus, there are a large number of unanswered 
questions even in the seemingly simple problem of the dynamics, near threshold, of a crack 
front restricted to move in a plane. 

In this paper we study the dynamics of a crack front restricted to move in a plane, 
through a three dimensional solid with heterogeneities only in the local fracture toughness. 
The effects of both the long range interactions and the stress pulses are considered, and some 
of the questions raised above addressed. In the absence of the stress overshoots, we obtain, 
numerically, some of the exponents which characterize the transition from the stationary 
to the moving phase, check the scaling laws that have been predicted and compare the 
exponents with the analytical results obtained earlier by the 2 — e expansion [0]. We then 
extend the analysis to include the effects of sound travel time delays in the stress transfer. 
Finally, we treat the effects of the stress overshoots on the depinning transition. Both the 
dynamic stresses obtained from a scalar approximation to elasticity and sharp pulse-like 
overshoots are studied. 



A. Outline 

Before introducing the basic model and summarizing our main results, we give an outline 
of the paper. In section || the details of the models and the numerical methods employed 
are described. Section [III A| contains the results of the quasi-static model, where the stress 



transfer is instantaneous, while Section IIIB contains those in the case where there are 



acoustic time delays in the stress transfer. In Section the effects of various kinds of stress 
overshoots are explored. Finally the results and their possible implications are discussed in 
section [V|. The "no-passing rule" |Q| for these models, which plays an essential role in the 
analytical results, is discussed in Appendix Appendix ^ has the detailed forms of the 
kernels used in the numerical simulations. 
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B. Summary of Results 



The equation of motion for the crack front can be obtained by requiring energy conserva- 
tion at all points on the crack front. This implies that the elastic energy flux into the crack, 
which is a non-local functional in both space and time of the shape of the crack front as 
well as of the local velocity of the crack front, must be equal to the surface energy required 
to create the new crack surface, i.e., the local fracture toughness. The general linearized 
equation of motion for a crack front moving along the positive x direction has the form, 

dtfiz,t)= J dt'vJdz'Jiz-z',t-t')d,>f{z',t')-^[fiz,t),z]+£ (2) 
t'<t 

where z is the co-ordinate along the crack front, V denotes the principal part of the integral, 
f{z, t) is the deviation of the crack front from a straight one, 7 is a random variable associated 
with the random position dependent fracture toughness in the solid and S represents the 
driving "force" due to the applied load, G°°. The kernel J is non-local both in space and 
time. This non-locality arises from the long range elastic interactions and the sound waves 
which are emitted as the crack moves. Note that because the basic processes near threshold 
consist of sections of the front moving ahead and stopping — i.e., roughly step functions 
in time — we have chosen to write the stress transfer in terms of d^' f{z' ,t'), so that these 
jumps are approximately delta functions in t'. 

We will classify the models based on whether or not the kernel, J, is monotonic in time 
at every spatial co-ordinate, z. Monotonicity of the stress transfer plays a crucial role. It 
means that as a segment of the crack moves forward, the stress at all other points increases 
monotonically in time. This convexity property yields stringent constraints on the behaviour 
as shown in Appendix ^. It implies that a configuration of the crack which is behind another 
configuration at one time will remain behind the other configuration at all later times. This 
immediately leads to the conclusion that there is a unique critical load, C}.^ for monotonic 
models. 
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1. Quastatic Approximation 



We first consider the quasistatic approximation in which sound waves are neglected and 
the stress transfer is instantaneous so that the kernel is naturally monotonic. In this case, 
the basic phenomenology is well known [0,0. As the load, G°°, is gradually increased, 
segments of the crack front will overcome the local toughness and jump forwards, perhaps 
causing other segments to jump and thereby triggering an avalanche which will eventually 
be stopped by tougher regions. We find that, similar to driven interfaces, etc. the 
avalanches show a power law size distribution up to a characteristic length, with larger 
avalanches being much rarer. The distribution of avalanche size — roughly the extent along 
the crack front of an avalanche — has the form 

Prob(size of avalanche > /) ~ ^P(V^-) (3) 

The cutoff length, defines the correlation length below threshold. As the threshold load, 
Gf^, is approached, the correlation length, ^_ diverges as 

e_ ~ (Gf -6-°°)-^-. (4) 

At the threshold, there is no characteristic length scale and the distribution of the avalanche 
sizes is a pure power law. From Eq.(^ and scaling relations between the exponents, we 
expect that the cumulative probability of the size of an avalanche being greater than /, as 
the load is swept slowly from zero to the critical load scales as 

^c,qs 

J rfG-^p(//e-) ~ y; (5) 



this is in agreement with the numerics within error bars. 

As the load increases above the critical load, the crack front begins to move with a mean 
velocity, which the monotonicity implies is unique. The velocity scales as in Eq.(|I]), with 
the velocity exponent, determined from our numerical simulations, 

/? = 0.68 ± 0.06. (6) 
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All quoted error bars here and henceforth are one-cx errorbars from fits. Just above G^** the 
motion of the front is very jerky with fluctuations in the velocity correlated up to a distance 
which diverges as one approaches the threshold from above as ~ {G°° — 
The exponents and z/_ will be equal i.e., = z/_ = u, if there is only one divergent 
length scale in the problem, as predicted by the renormalization (RG) analysis Assuming 
this two sided scaling, we can obtain the correlation length exponent, via finite size scaling, 
from the dependence of the variance of the critical load on the size of the system as 

1/ = 1.52 ±0.02. (7) 

At threshold, the crack front is self-affine with correlations 

([/(z,t)-/(z + r,t)f)~r2f (8) 

where () denotes the average over the randomness. The roughness exponent ( is found 
numerically to be 

C = 0.34 ±0.02, (9) 

in excellent agreement with the 2—e expansion prediction ( ^ 1/3. 

The dynamic exponent is found from the duration of avalanches as a function of size /; 
they typically last for 

Tl ~ I' (10) 

with 

= 0.74 ±0.03. (11) 

The exponent identities predicted from the scaling and RG analysis 

^={z- Qv (12) 
^ (13) 



1-C 

are found to be satisfied, so that there are only two independent exponents, say C, and z, 
characterizing the transition. 
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2. Effects of Sound Travel Time. 



For a model with a monotonic kernel but with the stress transfer delayed by the sound 
travel time, we argue that the static exponents u and ( are identical to the corresponding 
quasistatic case. Also, for every manifestation of the randomness, the critical load for this 
model, G^f, is exactly equal to that in the quasistatic approximation, G^^. However, the 
dynamic exponent for this model is predicted to be z = 1 exactly. Since the exponent 
identities Eq.(|T^ and Eq.(|T3|) also hold for this model, we obtain (3=1, which is consistent 
with the numerical results. 

3. Sound Waves and Stress Overshoots. 

The inclusion of the effects of sound waves leads to non-monotonic kernels. These result 
in the stress at points on an advancing crack front being non-monotonic in time, which 
substantially changes the physics. We have considered two types of non-monotonic kernels, 
the first arises from a scalar approximation to elasticity theory and the second is a simpler 
one characterized by sharp pulses superimposed on the time-delayed stress transfer. In both 
cases we find that the overshoots in the stress are relevant at the depinning transition, 
causing large avalanches to run away and changing the nature of the transition from the 
pinned to the moving phase. 

The model with sharp pulses involves non-monotonic kernels of the form 

Jsp{z, t- a, 7) = e(t - \z\)/tiz^ + a5{t - \z\)/\z[< (14) 

For a = 0, there are no stress pulses and the model reduces to the sound travel-time delayed 
model and hence Gc{ol = 0,7) is identical to the threshold force for the quasistatic model, 
G^*. We find, both from analytic arguments and from the numerics, that for small positive 
a and fixed 7 > 1/2, the threshold load, Gc{(y,'j), decreases with increasing a as 

{Gr-GT{a,j))r^a\ (15) 
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This behaviour is controlled by the relevant eigenvalue for the overshoot perturbation at the 
quasistatic depinning fixed point. 

In a scalar approximation to elasticity theory, the stress overshoots have long tails in 
time. In addition, the rough crack front will affect the propagation of the stress pulses due 
to the non-linearities neglected in Eq.(0). We argue that the basic features found in the 
sharp stress pulse model still obtain, in particular that the stress overshoots are relevant 
and change the nature of the transition. Numerical results using an appropriate class of 
kernels support this conclusion. 

For real elastodynamics appropriate to a tensile crack, the stress transfer kernel, for a 
fixed z, is found to be initially negative, when the longitudinal sound waves arrive, and then 
change sign when the Rayleigh waves arrive. The stress peaks due to the Raleigh waves 
are similar to those in the scalar elastic approximation and we believe that they will have 
similar effects in decreasing the critical load. However, the more complicated nature of the 
stress transfer suggests that the depinning transition of tensile cracks may involve essential 
additional physics. Some tentative ideas in this direction are discussed at the end of the 
paper. 

With or without the additional complications of the full elastodynamic stress transfer, the 
nature of the transition between a static and a moving crack front in the presence of stress 
overshoots is not resolved by our numerical or analytical results. The simplest scenario, 
which appears to be supported by the numerics, is a "first-order" transition with hysteresis 
from a pinned phase to a state with non-zero velocity. This may well be the correct scenario, 
but possible concerns and other possibilities are discussed in section V. 

II. MODELS 

In this section we discuss the equation of motion for a real tensile crack and various 
approximations to it that we will study. 
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A. Geometry and Equation of Motion 



We denote the plane in which the crack is confined y = 0, with the crack open in the 
region x < F{z, t). We assume that F{z, t) is a single valued function of z so that, the curve 
X = F{z,t) describes the location of the crack front. Since the crack is planar, the fracture 
surfaces that it leaves behind are, of course, smooth. The geometry is shown in Fig. 

The vectorial displacement field u, satisfies the equations of elastodynamics 

pd^Ui = djaij (16) 

with (Tjj the stress tensor. The displacement field, u{x,y = 0^,z), has a discontinuity 
across the crack surface while the normal stresses, criy{y = 0"*^), must vanish on the crack 
surface. For a crack with purely tensile loading, only Uy will be discontinuous and will have 
a \Jf{z, t) — X singularity at the crack front with an amplitude proportional to the local 
mode I stress intensity factor Kj{z,t) |^. As long as the crack remains planar, symmetry 
under y —y implies that the loading is purely mode I, so that we will simply use K = Ki 
rT| . We consider the system under a static load applied far away so that for a straight crack 
at rest (i.e., F{z,t) = const), K = K°° = const. 

As the crack front advances, F — > F + 6F, an energy per unit area of the new crack 
surfaces exposed, r[x = F{z,t),z], must be provided to the crack front in order to fracture 
the solid; in an ideal quasi-equilibrium situation this is just twice the solid-vacuum interfacial 
energy density, more generally it is the local fracture toughness that includes the effects of 
small scale physics for which linear continuum elasticity is not valid. The fracture energy 
will be provided to the crack-front by a fiux of stored elastic energy per unit area of the new 
crack surface, Q{z,t, {F}), which in general depends on the past history of the whole crack 
front as well as its instantaneous local velocity ^ . The equation of motion of the crack 
front is obtained by requiring that the elastic energy released be equal to the surface energy 
required for fracture, i.e., 

g[z,t,{F{t < t)}] = T[x = F{z,t),z] (17) 
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for all z and t. The available energy , has the general form 

Q = A\v^{z,t)\G\z,t{F{t<t)}\ (18) 

where u_l is the local velocity normal to the crack front and G, which is independent of 
^{z,t), is the elastic energy that would be released at {z,t) if the crack had advanced 



adiabatically at that point, i.e., with ^{z,t) = JTO 
For a straight stationary crack, F = const, 

l-z/2 



G = = ^——{K°°y, (19) 
E 



with E the Young's Modulus and u the Poissons ratio [Q. When the crack advances at a 
non-zero velocity, not all of the released elastic energy is available for fracture; some fraction 
of it goes into the kinetic energy of the moving material very close to the front. The fraction 
of G available for fracture A[f_L(z, t)] depends only on the local normal velocity; it decreases 
from unity for small Vj_ and goes to zero for v± = c, the Raleigh wave velocity. For a straight 
crack in a system with uniform toughness, F, this leads to a monotonic 

v{G^) = A-i(G°°/F) (20) 

for greater than the Griffith threshold, i.e. G°° > F. When G°° is smaller than the 
Griffith threshold we assume that the crack does not move, i.e., once the solid breaks, the 
crack does not to reheal (This is in fact observed in most situations, with the absence of 
rehealing due to plastic and other irreversible deformations at the crack tip). The velocity 
of the crack is thus constrained to be positive. [jl3[ 

We are interested in the behaviour near the depinning transition at which the crack 
starts to advance. We will use Eq.(|l^) as the starting point of our analysis of the dynamics 
of the crack front at this transition. The fracture toughness, F, in a heterogeneous solid , is 
a position dependent quantity, which we write as 

T{x,z) = Fo[l + 'yix^z)] 
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with Tq the mean value of the fracture toughness and ro7(x, z) the variable part of the 
fracture toughness which we will take to have a zero mean and covariance given by 



{'j{x, z)'j{x' , z')) = T{x — x', z — z') 



(21) 



with a function T which is, generally, short-ranged in space. 

The available energy, is a complicated nonlinear functional of the crack shape. In 
order to make progress, we will expand the position of the crack front in powers of the 
deviation, f[z, t), away from a straight crack. The position of the crack- front can be written 

as 



Fiz,t)=Fj + fiz,t), 



(22) 



where Fj is the original length of the crack which is assumed to be much larger than the 
scales of motion of the crack front so that the applied stress intensity factor, K°°, does not 
increase significantly as the crack advances. Thus, the stored elastic energy available to the 
crack front can thus be written in the form 

g = G°^[l + g{z,t,{f})] 

where G°° is that for a straight crack of length Fj for the given external load. If df/dz is 
small, so will be g. To linear order in /, g can be written as, 

g = -P0f 



where P is a kernel and ® represents a convolution in space and time. For a tensile crack. 



the Fourier transform of P is [111 [14 



P{k, u) = { 2Jk^ - Luyc^ - Jk^ - Luya^ + 7^f dWI{W, k\ uj') } 



with 



-UJ 



In 



(23) 



(24) 
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the contour integral circling in the counter clock-wise direction the cut in the complex W 
plane that runs from W = b"^ to W = with a and b the longitudinal and transverse sound 
velocities respectively; the argument of the logarithm in Eq. ( |24D the function whose zero, 
Wq, determines the Raleigh wave velocity via Wo = c^; and u u + iO needed to define all 
cuts e.g., 



for cj^ > c^/c^. Thus, 



\/k'^ — io'^/c^ = —isign{ijj)\\J uj"^ / — k'^\ 



p{k = 0,cj) = -iuB 



(25) 



(26) 



where i? is a positive number. 

To linear order in /, the equation of motion can be written as 

with the constraint that the local velocity of the crack front be positive and 

-To 



(27) 



(28) 



which acts like the applied driving force on the crack front. From the general structure of 
the energy release, from Eq.(0), we can separate the kernel P into sum of two terms, 
one which just depends only on the local velocity of the crack front and the other which 
depends non-locally on the shape of the crack front at all prior times. Thus, we can express 
P as 



P{k,uj) = -iuoB + \k\P{uo/\k\) 



(29) 



where |/c|P(co'/|A;|) is the non-local part which vanishes as a; — >■ oo. The equation of motion 
can then be written in the form 

Bdtf=\v I Jiz-z',t-t')dt'fiz',t')-^{z,t)+£\x 

[ z',t'<t ] 



e 



V I Jiz- z',t-t')dt'fiz',t')-^{z,t)+£ 



z',t'<t 



(30) 
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where 

J{z,t) = J J e^'^^— *[-|ip(.;/|fc|)] (31) 

and V denotes the principal part of the z integral. The Heavyside step function, B, constrains 
the velocity of the crack front to be positive and will not be written out explicitly henceforth. 
The kernel J is readily evaluated from Eq.(p3|) to be 



ate{at-\z\) 2ctQ{ct-\z\) 
J{z,t) — ——77315 ZvTTP?. + 



7Tz^{aH^ - z^y/^ 7iz^{cH^ - z^y/^ 
1 r t \/ I W I W 



62 y V V 52 



where the branch cuts are defined as previously. 

The stress transfer kernel Eq.(|32D is rather complicated. It is plotted as a function 
of t for a fixed z in Fig0 for a Poissions ratio, u = 0.25. At the arrival time of the 
Raleigh waves, J diverges as l/(ct — z)^/^ and then decays slowly to its long time value, i.e., 
J{z,t 00) the static stress transfer kernel. Although the negative stress precursor 

to the stress peak that occurs for z/a < t < z/c may well be important, for our primary 
purposes here, we believe that the stress peak is the more important feature. It is therefore 
useful to study a somewhat simpler model which has a similar stress peak. 

We choose to study a scalar approximation to elasticity theory. In this approximation the 
displacement field in the solid is take to be a scalar field 0, satisfying the three dimensional 
scalar wave equation 

- VV = (33) 

The displacement field (p has a discontinuity across the crack surface while the normal 
derivative dyip{y = 0^) (the "stress") vanishes on the crack surface. We shall refer to 
this model as the scalar model. Under the external load, the displacement field has a 



F(z^t) — X singularity at the crack front proportional to the scalar stress intensity factor 
K{z, t) as for real elasticity. The corresponding kernel P can be written in Fourier space as 
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Pscalar = " ^V^' (34) 

and the stress transfer kernel is 

ctQ{ct-\z\) 

We see that the stress peak, the long time tail and the static stress transfer kernel J{z, t — *■ 
oo) are all of similar form to the real tensile crack case. 

From the kernels Eq.(32) and Eg .(1351) we see that for both the tensile crack and the 



scalar model, sound waves yield non-monotonic kernels which lead, in response to a jump 
of one segment of the crack front, to ephemeral overshoots of the stress above the eventual 
static value. The magnitude of the overshoots will be governed by microscopic factors such 
as the microscopic response time of the crack front and acoustic damping processes which 
can be incorporated by the replacement of uj"^ by Vt^ = — in Eq.(El) or Eq.(|3^) with 



an acoustic relaxation time. We are interested in how these overshoots affect the dynamics 
of the crack front near threshold. But in the limit that all crack disturbances move along the 
crack slowly compared to c, we can neglect the effects of sound waves, and the transfer of 
stress will be effectively instantaneous yielding the quasistatic model with the kernel given 
by i 

Jqs(^,t) = ^e(t). (36) 

This quasistatic model we study first, its possible regimes of validity are discussed in Section 
V. 

In order to separate the effects of sound travel time delays from those of stress pulses, 
we also consider a model with monotonic stress transfer characterized by the kernel 

J,A{z,t) = ^Q{t-\z\). (37) 

This kernel is similar to the quasistatic kernel except that the stress transfer is not instan- 
taneous. Finally, in order to separate the effects of the maximum of the stress peaks from 
those of their tails, and to make the analysis of the stress peaks more tractable, we study a 

16 



kernel with sharp pulses defined in Eg. ([M]) . In both of these artificial models, the velocity 
of signal propagation has been set equal to one. Snapshots of the stress pulses when the 
crack front at z = is moved ahead by a small amount at t = and held there, are shown 
for the various models in Fig0. These are just plots of the respective J{z,t) for a fixed t. 

In general, the stress transfer along the crack front will depend on the shape of the front 
due to non-linear terms in the expansion of ^[{/}] in powers of /. Throughout this paper 
we will ignore these. We can justify this approximation for the quasistatic case for which we 
have an analytic understanding, and believe that it should generally be valid on long length 
and time scales as long as the crack front roughness exponent C < 1; i.e., that the crack 



front looks straight on asymptotically long scales. |11 



B. Numerical Implementation 

We are interested in the behaviour of the crack front near to when it begins to move. 
Below and just above the critical load, the motion of the crack front is very jerky and 
segments of the crack front move ahead and then get stuck in a tougher region. The basic 
minimum length, time and increment in the crack front position scales of these processes are 
set by the length scales of the toughness variations, the coefficient B in Eq.(^) etc. Thus 
to understand this behaviour, it is natural to simulate the crack front in a manner which is 
discrete in space, time and position. 

Our simulations were done on a lattice which is periodic in the direction, x, of crack 
advance, but for each z = 0...L — 1, the co-ordinate along the crack front, the row of points 
are shifted by an independent random amount, b{z) with < b{z) < 1. The allowed values 
of the crack front position are thus 

f{z) = b{z) + n{z) (38) 

with n{z) integers. This avoids the possibility of phase locked behaviour in which the crack 
advances in a relatively uniform manner characterized by all points advancing by one before 
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any advance again. Periodic boundary conditions are used in the 2;— direction. We have 
chosen the sound velocity in our models to be unity corresponding to one lattice spacing 
along the crack front per time step (all the models we study numerically have only one sound 
velocity). 

At each lattice point an independent value of the random fracture toughness, 'j{x,z), 
is picked from the interval [0, 1.5]. This range is chosen so that the variations in 7 are 
comparable to the force on points of the crack front on each other which are 

9{^, t) = E E ^(1 1^ - ^1 M - t') [/(^'' i') - /(^'' - 1)] (39) 

z'=Ot'<t 

where, with the periodic boundary condition in z on a crack of length L; 

1 1 2; — 2;' 1 1 = min{\z — z'\,\L — \z — z'W), (40) 

is the shortest distance between z and z'. The stress transfer kernels J are modifications of 
the continuum kernels of interest with the stress pulses designed to die away smoothly after 
going through the system once. Thus, although there is a long ranged history dependence 
in all but the quasistatic model, we need keep track of the history of the interface only up to 
a time corresponding to the sound travel time through half the system, i.e., a time of L/2, 
where L is the system size. Thus, 

J{z, t > L/2) = J{z, 00) = Jgs{z) = j^. (41) 

for z ^ 0. The sum in Eg . (p9D over —00 < t' < t — L/2 — 1 can thus be replaced by 
J2 J{\\z ~ z'\\)f{z'jt ~ L/2 — 1). Care must also be taken with the "self-interaction" piece 

z' 

z = z' which represents the "principal part" in Eg. (|30|) . In the guasistatic case (and generally 
for long time) 

j(^ = o,t) = -E4(lkll) (42) 

so that if the crack moves uniformly there are no changes in g. More generally, in particular 
for the artificial models, J{z = 0,t) involves some arbitrariness. To preserve the mono- 
tonicity of the sound-travel-time delayed model, a certain choice is reguired. This, and the 
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detailed form of the various J(z, t) used are specified in Appendix ^ The evaluation of 
the elastic force at each time step is done in Fourier space using the FFT algorithm and 
hence the time for computation at each time step of the evolution of the interface scales like 
L^logL. 

The "driving force", S, is the forcing parameter and as this is increased the crack begins 
to move. When the total "force" at a point z on the crack front is greater than the random 
part of the fracture toughness there, the crack front at z advances by one lattice constant, 
i.e., 

f{z, t+l) = f{z, t) + e[giz, t) - 7(^, t) + £] (43) 

where the lattice constant is set equal to one and 6 is the step function. 

These discrete automaton models for the crack front are expected to capture the physics 
at threshold of the corresponding continuum models at length scales long compared to the 
correlation length of the random toughness. Direct evidence for universality is provided by 
extensive numerical simulations on charge density wave models [|T5| in two dimensions which 
have found universal behaviour for smooth and piecewise continuous pinning forces as well 
as for discrete cellular automata analogous to the one defined above. We expect the same 
to hold here. 

III. MONOTONIC MODELS 

A. Quasistatic Model 

We first consider the quasistatic model. In this approximation, the stress transfer is 
instantaneous and the linearized continuum equation of motion of the crackfront takes the 
form 

dj{z,t) = -V I dz'l^^^j^^^^p^-j[f{z,t),z]+g (44) 
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Before presenting the numerical results from which we determine the values of various 
critical exponents characterizing the depinning transition, we give, following references 
scaling arguments for several identities between the exponents and bounds on them. 



1. Exponent Identities and Bounds 



In the moving phase, the crackfront will be reasonably smooth at scales larger than the 
correlation length ^, indeed at large scales / ~ f t and hence the random toughness, 7(/, z), 
will act essentially like white noise and hence {{f{z, t) — /(O, t))^) ~ ln(z). On scales smaller 
than ^, the front will be rough with t) — /(0,t)| ~ \z\'^ . 

In order for the crack motion to be smooth on larger scales, each segment of the crack 
of length ^ must take about the same time r ~ the correlation time, to move through 
each distance In the region which a segment of length ^ passes through in time r, there 
are ^^^^ random values of the local toughness. This means that the force per unit length 
needed to pull the crack segment through this region must vary from region to region by at 



least of order the random variation in the toughness averaged over this region, i.e., 1/^^^+^ 
by the central limit theorem. Thus the course grained toughness variations at the scale ^ 
are 



The force on the segment from the external load and the rest of the crack must be just strong 
enough to overcome these random variations. If these forces were too strong, the segment 
would move more smoothly implying that it must have been longer than ^ by definition. 
On the other hand, if they were too weak, the segment would not move at all in some time 
intervals of length r and thus it must have been smaller than ^. Since the mean external 
load at which the segment moves is 




(45) 



G. 



c — 



(46) 



this implies that either 



20 



- Ge ~ > (47) 



I.e., 



- > (48) 

or the force per unit length from the neighboring sections of the crack, Gn, are comparable 
to 6T^. The latter is dominated by nearby segments so that 

~ 1 dz^ ~ (49) 

yielding 
and hence, 

C > 1/3. (51) 

If there is only one basic scale of the forces near threshold, as simple scaling would suggest, 
then we should expect that 

- ~ G„ ~ ST^. (52) 

The rough equality of the typical force per unit length, Gm, of a segment of length ^ on a 
segment a distance ^ away and G°° — Gc thereby yields the scaling relation 

1 



(53) 



1-C 

Similar argument can be used below threshold implying that the correlation length exponents 
on the two sides of the transition are equal We will derive the relation Eq.(p3D more 
directly below. 

The bound on z/. Eg . (|48|) , comes from an argument similar to that by Harris |jl6| for 



equilibrium phase transitions and established more generally in reference fl^. We can also 
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derive it by considering some segment of the crack front of size ^_ below the threshold, 
loosely defining = Gc + to be the local critical force at which this segment begins 
to move. Below threshold, segments of size of order ^ have a substantial chance both of 
having already moved or not having moved yet while was increased to its present value. 
Thus the variations in must be comparable to Gc — G°° . This is similar to what was 
argued for above Gc, but here it does not rely on as many assumptions about scaling. Since 
5T^ > we obtain Eq. (|48|) just as from above threshold. 

We now obtain scaling relations for the distribution of avalanche sizes below the threshold 
loading. Following Ref. [Q] we conjecture that the distribution of avalanches as G°° is 
increased slightly has a scaling form 

Fraction of avalanches with size>/ when G°° G°° + dG°° ^ Ji^P^^I^-) i^^) 

at a given external load, with ^_ ~ (Gc — G°°)^'^. Following the same reference we obtain, 

K = l-l/y (55) 

from the increase in the mean position as Gc is approached. Now consider the probability 
distribution of the size of all avalanches that occur on sweeping the load from zero to the 
threshold load. 

Fraction of all avalanches with size >l ^ j -^p{l/i-)nA{G'^)dG°° ~ i-^-i/u ^ ^-i ^^q-^ 



where the last equality was obtained using the scaling relation Eq.(p5D and the observation 
that the rate of avalanche production, ?t,^(G°°), per increase in G°° goes to a constant at 
Gc. Thus, 

Prob(size of a given avalanche =/) ~ 1//^ (57) 

An exponent identity relating the velocity exponent, (3, to the other exponents follows 
directly from the picture discussed above of the moving segments of length ^. Since the time 
for a segment of length ^ to move ahead a distance which scales as C,'', is of order the 
velocity of the interface scales as We thus obtain 
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P = {z- Qv. (58) 

Another useful relation can be obtained by considering adding an additional "force", 
e{z,t), on the crack front. Denoting the resulting change 6f{z,t), we can define the polar- 
izability x 

6{f{k,u;)) 
^^^'") = 7^(^ 

in Fourier space. First consider applying a static force. In this case the additional force 
e(fc, uj) can be absorbed by redefining 

f(k,u;)^nk,u;)-'j^, (60) 

since the terms from the interaction of the crack front with itself in the equation of mo- 
tion will then exactly cancel the additional force e. The statistics of the random toughness 
variables, in this distorted frame, 

'y'{x,z) =-f{x + (f){z),z), (61) 

with (j){z) the Fourier transform of e{k)/\k\, will have the same statistics as the original ones. 
This is an important statistical symmetry of the system. It is the small angle form of the 
statistical rotational invariance. We thus have 5{f{k,uj)) = and hence 

X(fc,0) = ^ (62) 

exactly. On the other hand, on applying a low frequency spatially uniform force, e(a;), we 
should have 

^^u;x{0,u;) = ^^r^{G^-G,)^-\ (63) 
Generally, we expect that x{k, i^) will have the scaling form 

Xik, uj) ~ ^ -^^Xik^, (64) 

—tUJ 

with the form of the prefactors implied by Eq. (|5^ . In the static limit, X{kC,,u = 0) = 0, 
and 
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lim y(A:,c.) ~ - G^f-H^^^^^U^o- (65) 



Comparing Eg. (|62|) and Eq.(l65D we see that 

e(G'°° - G^f-^ 1 

But since ^ ~ - Gc)"'', we have 



(66) 



j3-l-v{z-l)=Q. (67) 

Using the expression for (3 from Eg . (^8]) we again obtain Eg . (|53D . As noted above, this 
simply relates the force of length ^ segments on each other to 6*°° — Gc- 

We thus have two independent exponents from which the others can be obtained. In 
addition, from Eg . (^81) and Eg .(^81) we obtain the bounds 

C > 1/3. (68) 

and hence, 

u > 3/2 (69) 



All of the exponent identities and the form of scaling functions such as Eg.(]6^) have been 
derived from a renormalization group expansion about two dimensions which is the critical 
dimension for the depinning transition of manifolds driven through random media with long 
range interactions decaying as l/r'^+^, — i.e., \k\ in Fourier space [l^. The analytical results 



from the d = 2 — e expansion are compared with our numerical results in the next subsection. 



2. Numerical Results 

We now present the numerical results for the guasistatic model from which we obtain 
the values of the various exponents. As discussed earlier, we simulate a discretized version 
of this eguation, here simply 
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f{z, t + l) = f{z, t) + - l{z, t) + £] (70) 

where L is the system size, to study the dynamics at threshold. We start with a pinned 
configuration of the crack front that is as close as possible to being straight and gradually 
increase the load until the most weakly pinned point becomes unstable and jumps. This 
point, in turn, may pull along other points on the crack front due to the elastic interactions, 
causing an avalanche. During the avalanche, the load is held fixed. Once the avalanche 
subsides, the load is increased once again until another point becomes unstable, and so on. 
A series of avalanches, as the load is gradually increased in this way is shown in Fig[Q. A 
space-time plot of one of the large avalanches is shown in Fig[^ ]. 

Defining the size of an avalanche is somewhat problematic. We have chosen to define 
it as the number of distinct points on the crack front that move during the course of the 
avalanche. Note that various other ways of defining the "size" of an avalanche along the crack 
front by e.g., its "moment of inertia" about its center of mass or by its maximum extent have 
problems because of the power law tail of the interactions which can trigger some jumps far 
away. In addition, periodic boundary conditions would complicate a definition. To study the 
statistics of many avalanches, the avalanche sizes are binned in powers of two. To measure 
the dynamic exponent z, statistics of avalanche sizes versus their durations are collected for 
all the avalanches that occur as the load is increased from zero to the critical load. The 
log- log plot of the number of avalanches in a given bin against the bin size is shown in Fig. ^ 
and a linear fit gives us a slope of 2. 14 ±0.3 in agreement with Eq.(p7D but with large errors. 

Figure[0 shows the plot of the mean duration of avalanches in a bin, Tun versus bin size, 
Ibini for a system of length 1024. From the slope of the log-log plot we determine 

z = 0.74 ±0.03. (71) 

As is generally true, one must be very careful not to take such statistical estimates of uncer- 
tainties in exponents at face value due to the existence of corrections to scaling. Fortunately, 
in our case, the 2—e expansion provides an estimate of the leading correction to scaling ex- 

25 



ponent; calculations of Narayan yield the leading irrelevant eigenvalue at the critical 
fixed point to be approximately — £:/3 ~ —1/3 in our case with e = 1. We thus fit the data 

to the form Tbin = ^-1/3 and find that Az << 1 and hence this fit gives the same value 

of z to within error bars. 

From the £— expansion []T|,0], it was found that 

z^l- 2£/9 + 0{£^) ^ 7/9 ^ 0.78 (72) 

for e = 2 — d = 1. U we neglect the O(e^) and higher terms, this agrees with our numerical 
result within error bars. 

In a finite system, there is some ambiguity in the definition of the critical load. For 
example, if the system extends very far in the direction of motion, the whole crack front 
would typically move from its initial position but get stuck in a rare tough region far away. 
In order not to bias the results by choice of the system extent in the direction of motion, we 
define the critical load as the load at which every point but one on the interface has moved 
at least once. For a large system we find 

G, ^ 0.97. (73) 

Note that the random forces, the critical driving force and the nearest neighbor elastic forces 
are very comparable. 

Right at threshold the crack front is found to be self-affine as expected. Figure [H shows 
the plot of the power spectrum of the crack front, {\f{k)\'^), at threshold, as a function of the 
wavevector for various system sizes ranging from 4 to 4096. We expect the power spectrum 
to be for small k. The best fit to a straight line is shown in Fig|^ for a system size 

of 4096 averaged over 1000 samples. The slope of this line gives us 2^ + 1, from which we 
determine 

C = 0.34 ±0.02. (74) 

Surprisingly, even at very small wavevectors, the power spectrum still looks linear on a log- 
log plot and we do not see significant finite size effects even at the wavevector corresponding 
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to half the system size. Since the data in Fig[^ have very small statistical uncertainties, 
we can try to fit the power spectrum to the form Ck~^'^'^^^\l + A^A;-"^/^) using corrections 
to scaling. The coefficient turns out to be very small and we obtain the same roughness 
exponent with comparable error bars. This gives us some confidence in the estimate Eq.(^). 

Our result for ( satisfies — and may saturate — the bound C > 1/3- The prediction of ( 
from the e expansion is 

C = ^/3+o(£") (75) 

for all n i.e., there appears to be no corrections to all orders in e although "non- 
perturbative" corrections cannot be ruled out. Nevertheless, the bound Eg. (|68|) , the e ex- 
pansion result Eg. ([75|) , and the numerics suggest that perhaps ( may be exactly 1/3 although 
at this point we know of no solid argument that yields 1/3 as an upper (to complement the 
lower) bound. 

The correlation exponent, z/, can be obtained via the finite size scaling hypothesis by 
measuring the variance of the threshold load, Gc{L), as a function of the system size. As- 
suming that there is only one important length scale ^, the variance of the threshold force, 
{AGc{L))'^ scales with the system length, L, as 

iAGc{L)f ~ L'^/". (76) 

A direct fit to L-^/" of the plot in FigQ of the variance of the threshold load versus 
the system size, for system lengths ranging from 4 to 8192, leads to u = 1.80 ± 0.05, while 
a fit for system lengths ranging from 256 to 8192, gives u = 1.72 ± 0.12. But a systematic 
curvature can be seen. In light of the knowledge of the corrections to scaling, we can do 
better by fitting to the form l-xi-i from which obtain 

z/ = 1.52 ±0.02. (77) 

This fit is shown in Figure |n[] for systems of length 4 to 8192 lattice constants. In this case. 



as suggested by the data. Ay is not small and the fit indeed yields = 0.87. Note that 
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the expected scaling equality Eq. (|53D is obeyed, but only when the corrections to scaling are 
included. 

As the load is increased above threshold, the crack front begins to move. As for the 
critical load, we must be careful how we define the velocity for finite length crack fronts. If 
we choose a system of great extent, W, in the direction of motion, the front will tend to get 
stuck in anomalously tough regions; this effect will be more pronounced for small L. But 
since we are interested in the critical behaviour and we have a good handle on the scaling 
of / with L, we can instead choose systems of extent W ~ CwL"^ with periodic boundary 
conditions in the direction of motion. For monotonic models, the convexity then implies 
convergence to a unique steady state §]. Since at threshold A/ ~ L'' with a coefficient 
roughly of order unity, we choose Cw = 4. 

There is a complication that must be considered: due to the possibility of a pinned 
configuration for loads above that at which the last point became depinned, the minimum 
at which v > will sometimes be greater than our definition of Gc by a random amount 
whose distribution depends on Cw From scaling we expect 



GZr^ - ~ (78) 

By scaling, there will thus be a typical minimum velocity 

v^,^ ~ ~ (79) 

Note that the minimum velocity due to the discreteness of time is much less than this and 
hence negligible for large L. 



Figure |12] shows a plot of the mean velocity of the front as a function of the loading for 



a system of length 1024 from which we determine 

p = 0.68 ± 0.06. (80) 

The fit using the corrections to scaling leads to the same value of f3 within the error bars. 
Surprisingly, there do not seem to be substantial deviations for G > G^ 



^imovmg 
'^min 
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The £— expansion prediction is /5 ~ 2/3 if we use ( = 1/3 and z = 7/9; again there is 
reasonably good agreement between our numerical results and the e expansion although our 
error bars are larger for j3 than for z, v or C,. 

In our numerical results, the critical force and all of the coefficients in scaling laws 
are of order unity suggesting that we have no intermediate length, displacement or time 
scales and thus that the scaling should work well for all but small size samples. If we 
had chosen too narrow a distribution of random toughness 7(x, z), there would have been 
an intermediate length scale and this would no longer have been the case. Note that the 
renormalization group methods can be used to show that non-linearities associated with 
higher order terms in the expansion of ^[{/}] are, for quasistatic stress transfer, irrelevant 
for the critical behaviour. One could have guessed this since, from the homogeneity of ^[{/}] 
higher powers of / have an equivalent number of powers of gradients, so that C < 1 implies 
that they are irrelevant. 

B. Monotonic model with Time Delayed Interaction. 

In the previous section, we saw that the quasistatic approximation to stress transfer 
gave rise to a critical depinning transition with a dynamic exponent z < 1. This means 
that for large enough avalanches which typically occur only if the load is close enough to 
Gci the effective disturbance velocity of an avalanche of size I will be l/l^ times the basic 
microscopic velocity scale of disturbances set by the dissipative coefficient B in Eq.(^). 
Thus for sufficiently large /, the quasistatic avalanches will progress faster than the sound 
speed. This is clearly unphysical and in this regime, the sound travel-time delays in the 
stress transfer must play a role. In order to understand the effects of these and of stress 
overshoots separately, we study a monotonic model — i.e. with no stress overshoots — but 
with sound travel-time delays. The simplest form of this is to simply replace the Q{z) in the 
quasistatic stress transfer with Q{t — \z\). On a lattice the stress transfer kernel becomes: 
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ll^ll ll^ll 

Some care is needed in choosing the second, local, term in Eg. (pT]) . The natural choice 
would be to fix J{z = 0, t) by the condition that for all times, X) J{z, t) = 0. This condition 

z 

is satisfied for the quasistatic model and for the scalar model as well as for real elasticity. 
It ensures that for a straight crack, the instantaneous crack front velocity is a function 
solely of the instantaneous external load and independent of the past history of the crack, 
since for a straight crack the effect of the crack front interactions vanish at all times. If we 
made =0 here also, however, the model Eq.([8T|) would no longer be monotonic. 

z 

Rather, the stress at a point z, after a jump at the same point, would decrease in time as 
the integrated stress transferred to the rest of the crack increases; this would act like a stress 
overshoot. For now we will, therefore, give up the independence of a straight crack on its 
past history to preserve the monotonicity condition. 

We see that in this model the force on any given point on the crack front at any given 
time is always less than or equal to the equivalent force for the same configuration in the 
quasistatic model. Note also that for a crack front which is stationary after some time t, 
in the discrete-time periodic boundary-condition version of this model (see Appendix 
the force at all points on the crack front, will reach the quasistatic value by time t + L/2. 
As shown in Appendix |^, these conditions imply that if the external load is increased 
adiabatically from the same initial conditions, the time-delayed model, Eq. (|8l|) , has exactly 
the same static properties as the quasistatic model. Indeed, for a given realization of the 
random toughness, every finite avalanche that occurs in the two models will be identical, 
except for the times at which points on the crack front jump. Therefore, both the models 
will have identical threshold forces, and the exponents u and ( are then obviously identical 
to their quasistatic values. 

The space-time plot of a particular avalanche as one site is triggered by increasing 
the load, is shown in Fig. [|13| for the quasistatic model and of the identical avalanche in 



the sound-travel time delayed model starting from the same initial configuration of the 
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crackfront in a system of length 512. It is evident that the dynamics of the avalanches in 
the two models is very different 

The Fourier transform of the continuum version of the time-delayed kernel in Eq.(|8T|) is 
given by 

1 i i 11 

J{k, uj) = { {k + uj)ln\k + uj\ H — {k — uj)ln\k — uj\ + \k — u!\} + U (u) 

iuj n TT 2 2 

(82) 

where there is an ambiguity in the uniform k = part of the Fourier transform, [/{uj). It is 
clear that the dynamic exponent must he z > 1. Let us assume that z > 1 or more precisely 
that the characteristic time t ^ ^ near threshold. We are thus interested in the behaviour 
in the scaling limit in which |u;| <^ \k\. In this limit 

J(k,uj) ^-\n\k\+^-^ (83) 

TT iuj 

From the equation of motion, the response function defined in Eq.(^) is in the absence of 
the randomness 

1 1 

—iujB + iujJ{k,uj) — ia;^ ln( |i|) + |A;| 

for UJ k, since the ln|A;| term from J will dominate the u dependence. In mean field 
theory, this gives rise to times scaling with lengths as 

r~eine (85) 

In the absence of the In \k\, renormalization due to the random roughness would make B^s 
decrease with length scale in 2 — e dimensions. However, the ln|A;|, being singular, cannot 
renormalize. But it can feed into the renormalization of B. Following this through yields 

~ (86) 

in dimensions d = 2 — e with c an effective velocity of order e for small e, but presumably 
of order unity in our one-dimensional case. Thus we see that assuming z > 1 leads back to 
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z = 1 (87) 

which we beheve should be correct with no logarithmic corrections. Using the exponent 
identities Eg. (|58|) and Eq.(|53D we obtain 

(3=1. (88) 



A plot of the load versus the mean crack velocity is shown in Fig[0 for a system of size 64. 
It is very close to linear although the range is small enough that one cannot reliably extract 
f3. Note that because of the dependence on the past history, we are limited here to rather 
small samples. While the computations for each time step for the quasistatic model take 
a time of order LlnL, those for the sound-travel time delayed model take L^lnL and hence 
even for a system of length 64, the statistics are more difficult to obtain. 

IV. STRESS OVERSHOOTS: NON-MONOTONIC KERNELS 

We now turn to more realistic stress transfer kernels. We have seen, from Eq.[^ and 
Eq. , that bulk sound modes naturally lead to non-monotonic kernels of the stress transfer 



along the crack front. A regime of negative stress transfer, as occurs at intermediate times for 
tensile cracks, cannot, by itself, change the behaviour much from the time delayed monotonic 
models since, in the absence of stress overshoots, the static behaviour will again be identical 
to the quasistatic model. Thus the primary differences between the time-delayed model and 
more realistic models must be associated with the stress overshoots. The actual shape of 
the stress overshoot may be complicated by various factors including the effects of multiple 
scattering of sound waves off the crack front. Therefore, we would like to understand what 
features of the stress overshoot play a crucial role in the dynamics near threshold. To do 
this, we study simpler models and hope that the conclusions drawn from these models will 
help us understand the case of real elastodynamics. 
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A. Sharp Stress Pulses 



We first consider a simple model of the overshoots in which sharp stress pulse travels 
with the sound speed. We take the amplitude of the overshoot to decay as a power law of 
distance as it moves along the crack front. The continuum version of the kernel we study 
has the form 

Jsp{z, t; a, 7) = e{t - \z\)/z^ + a6{t - \z\)/\z\^ (89) 

This kernel reduces to the previous case of the monotonic time delayed interactions when 
a = 0. 



In Figure [[Tq] a large avalanche that occurs on triggering the most weakly pinned site is 
shown for both the monotonic time delayed kernel and for the kernel Jsp{z,t; a = 0.5,7 = 
1.5) starting from the same initial configuration of the crack front and the same configuration 
of random toughnesses. We see that in the presence of the overshoot many more sites are 
triggered than for the monotonic kernel. Our data show that even for small overshoots their 
effects build and cause sufficiently large avalanches to run away. This causes the crack front 
to de-pin and start moving at a threshold load which is less than the one for the quasistatic 
stress transfer. 

We find that for any value of 7 and any non-zero a, the threshold load, Gc{<y, 7), is lower 
than Gc(« = 0,7) = Gf^. A plot of {Gf^ — Gc{oi,'^)) as a function of is shown in Fig||T6[ 



for various values of 7, including 7 = 00, for which only the nearest neighbor of a jumped 
site feels the overshoot. The shift in threshold load from Gc was obtained by averaging over 
the same set of random samples with and without the overshoot; this greatly reduces the 
error bars. The value of Gc can be accurately determined from the quasistatic model where 
the code is much less numerically intensive. For all the values of 7 studied, 7 > 0.5, we finds 
results consistent with 

Gf -Ge(«,7) (90) 
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for small a, a form which will be derived below. Thus the overshoot appears to always be 
relevant at the quastatic depinning fixed point. 

One of the advantages of starting with monotonia models for which we have quite a 
detailed understanding, is that the effects of perturbations away from these can be analyzed 
using known scaling properties of the monotonic models. We would like to carry this out 
for weak stress pulses added to the time-delayed monotonic model; i.e., to consider the 
behaviour of the crack with stress transfer given by Eq.(^) for small a. In order to do this 
we first obtain the response to a single stress overshoot. 

We focus on a given space time point, at {Z,T) which we denote "A" and an avalanche 
that started at (0, 0) in the time delayed model without stress pulses. For simplicity we 
restrict consideration initially to 7 ^ oo so that only the nearest neighbor of a jump site 
will feel a stress pulse which will be a above the static stress. In order for A to be affected, 
one of the two neighboring sites of A, say Z — 1, must have jumped at time T — 1 producing a 
pulse; denote this space-time point "P" . Three conditions must be met for A to be affected 
by the stress pulse from Z — 1, i.e., for point Z to jump an extra time. First, A must be 
within a of jumping anyway for the stress pulse to have been able to trigger its jump. The 
force increment needed for individual sites to jump are uniformly distributed for sites near 
to jumping, so the probability of this is simply of order a. Second, the increase in stress at 
site Z after time T must be less than a, or else the site would have jumped again regardless; 
denote this condition "L"(for "later stress"). Third, the neighbor must jump at T — 1, i.e., 
the jump P must occur. 

The scaling properties of avalanches imply that within the space-time volume of a large 
avalanche, there are holes on all scales ( since fl << /^) and sub- avalanches of all sizes; see 



FigfTSla]. Consider a sub-avalanche that occurs in a region within a time r << T before 



P and within distance of order r from P that is mostly contained within the "backward 



sound cone" of P, as shown in Fig||T^. If P is not well enough "isolated" from the bulk of 
this sub- avalanche, then we must consider a smaller sub-avalanche until we find the largest 
r such that there is a sub-avalanche with size of order r from which P is "isolated" in 
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space-time by order r ( if this does not exist, then it is highly unhkely that condition L can 
be satisfied). By isolated, we mean that there are few (or no) jumps on or outside of the 



backwards sound cone of P — "P's cone" — within a time r. This is illustrated in Fig[]T7 
The maximal such sub-avalanche which we denote "5"' typically has of order r^^^ jumps of 
which a fraction 1/r, i.e. r^, are on P's sound cone. Since these are typically a distance 
^ r from Z, they each cause the stress at P to increase by order l/r^. Thus the isolated 
space-time point P will typically have a stress increase of l/r"^^^ ( between time T — 2 and 
T — 1) due to the sub- avalanche 5", the probabihty of a jump at P given S is of this same 
order. 

The condition that P is isolated from S is not a stringent one. Any jump in 5* outside of 
A's cone ( which is almost identical to P's cone) will cause a stress increase, A^, at Z after 
time T. Since there will be of order t^^^ such jumps, each yielding a stress increase at Z of 
order will be of order . This stress increase will be larger than a, and hence 

violate condition L, unless r > \/a^/^^^^'> (note that obtaining this condition by either a 
smaller sub-avalanche with anomalously few jumps outside of A's cone or by P being less 
isolated from S is very unlikely). The probability of both satisfying the condition L that the 
stress increase after T be small and having a jump at P, is thus controlled by the smallest 



~ (91) 

for which L is satisfied with reasonable probability. The relevant subavalanches S are thus 
of size of order Tq, and occur in a space-time region of extent x near P. Since L is then 
quite probable 

Prob(jump at P and L|size (S)~ r^) ~ l/r^-^ ~ a^+i/^^-^), (92) 

which is simply proportional to the increase in stress at P due to S. 

What is the chance that there is such an appropriate sub-avalanche S in the space- 
time region within of P? If there is any activity in this region, then it should include 
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subavalanches on all scales, but we also need this to be the last activity in this region (or 
else L will be violated). If T is of order the duration of the full avalanche, the last activity 
in this region could occur anywhere within a time of order T. Thus 

Prob(last activity being within of T) ~ Tq,/T. (93) 

We obtain the probability of a pulse-triggered extra jump at A by combining all the factors 
from Eq. (|92|) , Eq. (|93D and the probability a of site Z being close to jumping again, yielding 

1 Ta 

Prob (pulse triggering extra jump at A\ avalanche of size ~ T) ~ -^a 

Tot ^ T 

(94) 



a" 



T 

using the Tq, given by Eq.(0). Since such an event could occur over a range of times of 
order T and the number of sites in the original avalanche is of order T, the total number of 
primary extra triggers caused by the pulses is of order 

Ni ~ a^T (95) 

The spatial density of the primary extra triggers is small so that each of them will cause 
roughly independent secondary avalanches under the dynamics without further stress pulses. 
The probability of large avalanches falls off slowly with their size up to the correlation length 
^, which we assume is bigger than T so that the original avalanche was not exponentially 
unlikely. This means that the total number of jumps in these secondary avalanches will be 
dominated by the largest one which has size 

/max ~ {N,Y'^ (96) 

Then the total number of secondary jumps caused by the primary extra triggers is 

M2 ~ («2y)(i+^)/'^ (97) 

using K = 
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Each of these secondary jumps has the potential of triggering more jumps due to pulses. 
The number of such secondary triggers, N21 will be much less than A^^i unless M2 ~ Mi; the 
original number of jumps. Thus, for fixed T, small enough a will cause a small number of 
secondary jumps, fewer tertiary ones, even fewer quaternary ones, etc. But if 

M2 ~ Ml ~ ri+^ (98) 

i.e. T ~ a~'^^^^~^\ the process will run away. Thus we expect that the stress pulses will 
become important when 

e ~ ~ a-'"-'-^^ (99) 

Using 1/z/ = 1 — C , this corresponds to a reduction in the critical load proportional to a^, 
i.e., of exactly the form Eg. (pO|) that provided a good fit ot the numerical data. 

Longer range pulses can be considered by a generalization of the above argument. It 
is found that the dependence is preserved if and only if / dzJp{z) < 00, with Jp{z) the 
peak pulse height (as a function of time) at distance z. Our results for 7 > 1/2 agree well 
with the predicted dependence of the reduction in the critical load. The marginal case 
7 = 1/2 is similar numerically and we have not explored smaller 7. 

For non-monotonic models, the no-passing rule discussed earlier and in Appendix 0, 
does not apply. Thus, at least for finite systems and for some length of time, moving and 
stationary solutions can co-exist. Indeed for any non-zero a, we expect that for loads in 
the range Gc{(y,'y) < G < G^^, both stationary and moving solutions should co-exist and 
the selection between the two will be determined by the past history. This effect is seen in 



Fig|jT8[ in which the velocity of a crack front of length 64, averaged over a cycle, is shown 
as a function of the number of times it has passed through a sample of extent W = IQ with 
periodic boundary conditions in the direction of motion. In monotonic models, convexity 
implies that at long times we would measure the same velocity in every pass as there is a 
unique steady state [Q. But with stress pulses, we see that the velocity changes with the 
cycle, there is no unique moving solution, and after a number of cycles the crack front can 
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suddenly come to a complete halt, as in Fig[jT8|, thereby directly illustrating the co-existence 



of moving and stationary solutions. In finite systems with periodic boundary conditions in 
both directions as we have used, whether all moving states eventually stop or whether they 
can survive indefinitely is likely to depend both on the sample and on the load. The question 
of what happens for infinite systems, we return to in the last section. 

B. Scalar Elastic Approximation 

We finally consider kernels of the form 

Jse[z, t, t) = — — -— (100) 

7lZ^[[t + Toj^ — Z^Y'^ 

which is appropriate for the scalar approximation to elasticity. Unlike the sharp pulse 
models the stress overshoot has a long tail in time. We have cut off the singularity at the 
sound arrival time by a time tq which crudely represents the microscopic response time. For 
To ^ oo, this model becomes the monotonic model with the time delayed kernel. 

We measure the threshold load as a function of tq, for large tq, and find the reduction of 
the critical load I^Gf — Gc(to)) is consistent with 7q ^ clS shown on a log-log plot in Fig[jT9| 



we do not, however, have an analytical argument for the exponent 3/2. Not surprisingly, we 
again find that the non-monotonicity of the kernel is relevant but with a larger eigenvalue 
that would be expected from the peak pulse heights, presumably because of the long time 
for which the overshoots are substantial. 

C. Velocity and Hysteresis 

For both the sharp pulse and the scalar elastic models, the stress overshoots lead to 
velocity versus load curves that are both very noisy and hysteretic. For the scalar elastic 



case, results are shown for various tq in FigpO]. As the load is increased, the velocity appears 



to jump to a non-zero value which is a function of the overshoot's strength and then jump 
back to zero again on decreasing the load only at smaller load. Thus the stress overshoots 
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appear to lead to a first order like transition, where the crack front jumps directly to a finite 
velocity from the pinned phase. Finite sized systems show hysteretic behaviour as shown 
in Fig[Q. Several cautionary remarks are, however, in order. First, getting statistics on 
the sizes of hysteresis loops for these non-monotonic models is numerically intensive. As we 
saw in Fig|]TB[, due to the non-monotonic nature of the models the crack front can come 
to a complete halt after passing through the system several times. The fluctuations in the 
mean velocity per cycle of the crack front increases with the size of the stress overshoot but 
for a given magnitude of the overshoot, decreases as the system size increases. The abrupt 
stopping of the crack front naturally leads to a large scatter in the size of the hysteresis loops, 
particularly for smaller systems. Thus, even though the computation time per step increases 
with system size as LHnL, the scatter in data of the hysteresis loop sizes for small system 
sizes makes it difficult to study the system size dependence of the hysteresis loops. Second, 
it is not possible to ascertain whether the crack front will eventually stop or not unless the 
load is above the of the quasistatic model; i.e., the largest load at which a static 

solution exists. Third, it should be noted that the minimum velocities both on increasing 



and on decreasing the load in Fig[^ are not all that much bigger than the quasistatic fmin, 
even for tq = 1 for which Gc has decreased by almost a factor of two. Thus, overall, it 
is not clear at this point which of the effects that are apparent in the numerics for these 
non-monotonic models are finite-size effects and which are indicative of the behaviour of 
much larger systems. We return to this issue at the end of the paper. 

V. DISCUSSION 

In this last section we compare our results on the dynamics of planar crack fronts with 
other work and discuss various open questions. 
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A. Quasistatic Limit 



In the absence of sound waves, long range elasticity leads to a non-local but monotonic 
stress transfer kernel in the equation of motion of the front. The transition from the pinned 
to the moving phase is then second order and there is a unique moving solution above 
threshold. There are two independent critical exponents in this case and the numerical 
results we obtain using a discrete model are in good agreement with those from the e- 
expansion. Note, however that it was necessary to include the effects of corrections to 
scaling both to get reliable estimates of the exponents and to verify the scaling laws. The 
dynamic exponent z and the velocity exponent [3 are both found to be less than one as 
predicted by the e- expansion. The roughness exponent C, is very close to the lower bound 
of ^ which may well be exact, although we have no solid argument for this. 

There have recently been two other numerical studies on quasistatic crack models with 
the appropriate long range interactions. The first, by Schmittbuhl et.al |21|], obtains the 
same value of the roughness exponent C, as we do to within error bars. However they obtain 
a dynamic exponent which is greater than one. It is not clear from their paper as to how this 
result was obtained and at this point the discrepancy is not understood. The second paper, 
by Thomas and Paczuski |^, obtains C ~ |- Their system sizes are large and it is not at 
all clear why the results should be so different. One possibility is the dynamical "updating 
rules". Thomas and Paczuski's are different and more unphysical than ours. Rather than 
increasing the force adiabatically to depin the most weakly pinned point on the crack front, 
they depin the point on the crack front that is farthest behind the rest. Whether this or 
some other difference is the cause of the differences we leave as an unresolved question. 

In addition to these numerical studies, there is a very recent experiment pSf in which 
two halves of a block of plexiglass which have been roughened and then pressure welded 
together are broken apart. The crack is thus confined to the plane of the original weld. The 
crack front roughness is measured while it is advancing at a very slow mean speed. With 
a rather limited range of data, the authors obtain C, = 0.55 ± 0.05. This would appear 
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to be inconsistent with our quasistatic result, but a systematic curvature appears to be 
observable in the data and it may be popssible to fit it reasonably well with C = | a 
correction to scaling (a form with the same number of parameters as an unknown 
(). Another possibility, however, is that the experiments are not really in the quasistatic 
regime. 

These experiments force us to address an issue which we have so far avoided: what 
determines whether (or in what regime) a crack will behave quasistatically? The basic 
criterion is not directly related to the average velocity of the crack. Rather, it is the speed 
of propagation of disturbances along the crack front - in particular this speed relative to the 
Rayleigh wave speed, c - that is the essential determining feature. But what determines the 
speed of propagation of disturbances is rather subtle. In general, if the materials that make 
up the heterogeneous medium are themselves resonably close to ideal elastic solids, then 
there is no natural parameter which would make the speed of propagation of disturbances 
along the crack front much slower than c, even if the heterogeneities are weak. But if there 
is substantial plasticity, creep, or other dissipative effects on the scale of the heterogeneities, 
then even on these mesoscopic scales the equation of motion of the crack front is not given 
simply by the "propagator" Eq.(pUD with B its ideal value of order 1/c. If the system is 
velocity toughening due to these small scale effects, i.e., the effective fracture toughness on 
the scale of the heterogeneities increases with velocity - or equivalently that the velocity 
increases more slowly with load than in an ideal solid - then the linearized dynamics is given 
by Eq.(^) with a larger value of B. In the limit of very large B, the large scale behavior 
of the crack front will be well approximated by the quasistatic model except very close to 
the onset of crack propagation where the cumulative effects of the stress pulses caused by 
a large avalanche will still cause it to run away. Since the exponent z is not much smaller 
than one, however, the crossover to fully dynamical behavior will occur only very close to 
the critical load. 

Other non-linear effects - in particular those associated with the local depinning of a 
section of the crack front caused by the advance of other sections of the crack - will also 
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affect how quasistatic the behavior of the crack front will be. A careful study of several 
experimental systems - including investigating whether or not the motion appears in bursts 
of activity when the crack is moving at a slow average speed - would appear to be needed 
to help resolve this and other related issues. 

B. Elastodynamic Effects 

Close enough to the critical load the quasistatic approximation always breaks down. 
Since in this approximation the dynamic exponent z < 1, the effective propagation velocity 
of the disturbance associated with a quasistatic avalanche of size ^ diverges as thereby 
becoming of order the sound speed sufficiently close to the critical load no matter how small 
the "bare" velocity of small scale disturbances. Thus elastodynamic effects must alter the 
asymptotic critical behavior. 

In order to understand the effects of sound-travel time delays, we first considered a 
simplified causal model in which the dynamic stress transfer is still monotonic. In this 
model, the static exponents (i.e. ( and z/ ) were found to be the same as in the quasistatic 
approximation because of the monotonic character of the stress transfer. But the dynamic 
exponent, z, became equal to one, the minimum value consistent with causality. The scaling 
identities then imply that P — 1 which is in good agreement with our numerical results on 
this sound-travel-time-delayed model. 

But the actual dynamic stress transfer along a crack is more complicated. Indeed, proper 
inclusion of the dynamics of the medium necessarily leads to a non-monotonic kernel in the 
equation of motion of a crack front. In particular, the stress that arises at a point on the 
crack front due to another section of the crack moving forward, generically rises to a peak 
before decaying to its long time quasistatic value. 

We have examined the effects of these stress overshoots and find that they are it always 
relevant at the depinning transition. Specifically, for sufficiently large avalanches the effects 
of the overshoots build up enough to make such avalanches run away. This causes the crack 



42 



front to move - at least by some amount - at a load which is lower than the quasistatic 
critical load, i.e., for loads at which there are still stable static configurations. (This can 
occur because, in the presence of stress overshoots, the "no-passing" rule which prevented 
moving and stationary solutions coexisting in monotonic models is no longer valid .) 

If the stress overshoots are weak, their effects will not be important until very close to the 
quasistatic critical load and there will be a wide regime of validity of the quasistatic results. 
They will eventually break down only when the correlation length exceeds a crossover length 
scale which has an inverse power law dependence on the magnitude of the stress overshoots. 

What happens when a large avalanche runs away? We have explored this by numerical 
studies using simplifed stress transfer kernels which include both sound-travel-time delays 
and stress overshoots. In order to investigate hysteretic effects, we have used finite systems 
with periodic boundary conditions in the direction of motion with the extent in this direction 
proportional to the cube root of the length of the crack front - i.e, the scaling of the crack 
distortions at the quasistatic critical load. In the absence of stress overshoots, moving 
configurations converge to a unique periodic state. But, due to the non-monotonicity, this 
need not be the case once stress overshoots are included. 

Nevertheless, for a range of loads between the point Gmn at which an avalanche runs 
away and the quasistatic critical load G^*, we find that the crack front usually converges to 
a state which is periodic in time with a period which is some multiple of the time to pass 
through the system once. If the load is then decreased to below Grun, the resulting moving 
state coexists with static configurations which are stable to avalanche runaway under small 
increases in the load. At still lower loads the behavior tends to becomes chaotic, with, in 
at least some samples, the crack eventually coming to rest only after passing through the 
sample many times as illustrated in Fig||18||. 

The data we have collected thus suggests hysteretic behavior with coexisting moving and 
stable stationary regimes coexisting in some range of loads from Gmn down to some lower 
critical load, Gstop- However, our numerical results indicate that the widths of the hysteresis 
loops are quite a bit smaller than the difference between the critical load Gmn and the the 
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quasistatic critical load. We would of course like to know whether the hysteretic behavior 
persists in an infinite system. Unfortunately, the numerics are rather slow in the presence 
of sound-travel-time delays. Thus, obtaining statistics for the sizes of hysteresis loops is 
numerically intensive and our results are far from conclusive. 

C. Possible Scenarios near Threshold 

In this last section we consider various possible scenarios for the behavior of large systems 
in the presence of stress overshoots. 

The simplest scenario is suggested by our data: As the loading is increased slowly, the 
stationary crack jumps to a non-zero velocity at load Gmove, but when the load is decreased, 
the crack does not stop until a lower critical load Ggtop- At Gstop the velocity could either 
drop to zero discontinuously, presumably the result of an instability of the moving phase, 
or continuously (as occurs in a hysteretic underdamped Josephson junction). If the load is 
changed suddenly, this could cause a jump from one phase to the other, perhaps even in 
the range above Gmove but below in which static configurations still exist . Note that 
the obvious guess, suggested by our numerics, is that G^nove — Gmm the point at which 
avalanches run away and become much bigger than their size in the absence of the stress 
overshoots. However, it is not obvious that this has to be the case: One could imagine a 
scenario in which the runaway avalanches eventually stop but only after causing the crack 
front to become much rougher than it would be under the quasistatic avalanches. A true 
moving phase might then only exist above a higher load Gmove- This appears rather unlikely 
and seems difficult to reconcile with our numerical simulations even though the latter may 
have been biased by our choice of scaling of the length and the extent of the finite systems. 

A second scenario is that, in the limit of large system size, the hysteresis loops we found 
numerically disappear and the transition becomes "first order " with a discontinuous velocity 
versus load but no hysteresis if one waits a long enough time for the crack to settle down. 

Finally, and perhaps most interesting, is the possibility that the onset of crack motion 
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could still be critical with the velocity rising continuously at a critical load and some kind 
of diverging correlation lengths as the critical load is approached adiabatically from above 
and from below. This would represent a new universality class of depinning transitions. A 
variant of this, with the velocity discontinuous but the transition still critical in the sense 
of diverging correlation lengths, is also conceivable. 

Which of the above scenarios obtains may well depend on aspects of the physics that 
we have left out of our numerical studies and theoretical analysis. For example, multiple 
scattering of elastic waves from the crack front will cause different regions of the crack front 
to see stress pulses that depend on the shape of the crack front in their vicinity and that of the 
segment which has moved. However, we conjecture that the general role of stress overshoots 
should not be qualitatively changed by multiple scattering since the long wavelength sound 
waves will not be strongly affected by the crack roughness unless C = I. 

A potentially more important effect is a consequence of vectorial elastodynamics. In 
particular, for a tensile crack the behavior may be complicated by the fact that the initial 
stress pulse caused by a section of the crack jumping forwards is negative with the stress 
only becoming positive when the Rayleigh waves arrive. If one hypothesizes a hysteretic 
velocity-load curve, then in the hysteretic region such a stress transfer kernel can support 
the coexistence of moving and stationary zones of the crack front with the boundary between 
the zones moving at a velocity s that corresponds to the zero of the kernel P{k,uj) ( for 
a Poisson's ratio of |, s ^ 0.94c). In the moving phase, an anomalously tough region 
could thus cause a stopping "shock" to move along the crack front. This might result in 
a complicated - and very rough - moving state involving large scale stopping and starting. 
What roles such shocks might play in the onset of macroscopic crack motion we leave as an 
interesting avenue for future study. 

Finally, we have totally ignored all effects of non-planar crack deformations. These almost 
certainly play a major role in many experimental situations and may well be important 
whenever the crack is not confined to a pre-weakened plane. 
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APPENDIX A: 

In this appendix we discuss the no-passing rule in the context of the quasistatic model 
and show that the monotonic model with the time delayed interactions and the quasistatic 
models should have the same static exponents. 

The elastic force in the quasistatic model can be derived from an elastic potential defined 

as 

V{{f])= J dk\k\\f{k)\' (Al) 

where k is the wave vector corresponding to the z co-ordinate and f{k) is the Fourier trans- 
form of f{z). This potential is convex in /. Following Ref. ^ consider two configurations of 
the crack front, f^{z,t = 0) and f^{z,t = 0) with f^{z,t = 0) > f^{z,t = 0) V z. The no 
passing rule states that this inequality holds for all times. This can be seen by noticing that 
if the two configurations were to pass, would have to first approach at some point 
z. At this point, the random fracture toughness and the external driving force would be 
identical for both the conformations of the crack front. However, the elastic forces at z on 

would be less than or equal to that on f^at z, due to the convexity of the potential. This 
prevents the passing of bby and hence the inequality f^{z,t) > f^{z,t) is obeyed at 
all times. The no passing rule also implies a unique moving solution for the crack front in 
the quasistatic model, p. 

Now consider the monotonic model with the time delayed interaction. We see that 

Mz,t) < Jqs(^,t) W{z,t) (A2) 

where Jtd{z,t) and Jqs{z,t) are the kernels describing the elastic interactions in the two 
models. This inequality holds when we define the kernel with the sound-travel-time delayed 
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interactions at 2; = as in Eg. (pT]) . Thus we see that if we consider two crack fronts, fqs{z, t) 
obeying the quasistatic equation of motion and /td(^, t) obeying the monotonic time delayed 
equation of motion, with /qs(-2, 0) = /td(-2, 0) V z, then following the previous argument we 
see that fqs{z,t) > ftd{z,t) V {z,t). Also since, as t — 00, Jtd{z,t) — * Jqs{z) V (z), we see 
that, if the load is below threshold, at the end of the avalanche fqs{z) = ftd{z). Thus if we 
start with the same initial configuration, the final positions of the crack front at the end of 
each avalanche are identical. Thus the static properties are the same for both the models. 

By defining the kernels as in Eq. (|4^) , we make sure that both our scalar model as tq — 00 
and the sharp pulse model with a = behave like the monotonic model with the time delayed 
interaction. 



APPENDIX B: 

In this appendix we give the explicit forms of the interaction kernels of the various models 
that we studied numerically. 

For the quasistatic approximation, the discretized version of the interaction kernel. 
Eg . (p6l) , with periodic boundary conditions in the direction along the crack front is given, 
as in Eg.(|l]), by 

= ^ (Bl) 

for 2; 7^ and 

4(^ = 0) = -E4W> (B2) 

where 

I |z| I = mm(|z|, |L — |z| I). (B3) 

In the case of the monotonic model with the time delayed interactions, the discretized 
version of the interaction kernel Eg. (P7|) is given by 
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J,,{z,t) = j^Q{t-\\z\\) (B4) 



for z ^ 0. As noted in the text and in the previous appendix, there is an ambiguity in 
defining the kernel at 2; = 0, and in order to preserve the monotonicity properties of the 
kernel we define 

Jtd{z = Q,t) = -Y,J,,{z) (B5) 

The other models we consider in this paper have stress overshoots which decay as the 
stress pulse moves along the crack front. To include these effects with periodic boundary 
conditions, we design our kernels such that the stress overshoots disappear smoothly after 
running through the system once, after which time the kernel equals the quasistatic kernel, 
Jqs at all points in space. Thus, for the sharp pulse model defined by the kernel Eq. (p9D, 
we choose the discretized kernel to be 



< L/2) = ^e{t-\\z\\) + J^^ JI'^^K -^^ (B6) 



and 

Jsp{z,t>L/2) = -\- (B7) 



\z 



for 2; 7^ 0. There is again an ambiguity as to how one chooses the kernel at z = and we 
have defined it to be 

Jsp{z = 0,t) = -max[J2 Jsp{z, t), Jgs{z = 0)]. (B8) 

Finally, in the case of the scalar model, the discretized form of the kernel, Eq. (p5D , is 
chosen once again such that the overshoots vanish after they have run through the system 
once. Thus, 

le{z, t < L/2) = i , (B9) 

and 
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Jse{z,t>L/2) = j^ (BIO) 

\z\r 



for z 7^ 0. As for the sharp pulse model, we define 



Jse{z = 0,t) = -maX^ Jse{z,t), Jqs{z = 0)]. (Bll) 
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FIGURES 




FIG. 1. Schematic of a planar crack propagating through a heterogeneous medium. The crack 
front X = F{z, t) and the free crack surfaces, which are flat, are shown. The applied mode I (tensile) 
loading is also indicated. 
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FIG. 2. Stress transfer kernel for a tensile crack is shown as a function of t at fixed distance 
z = 1. The longitudinal sound velocity is set equal to one and the Poisson ratio chosen to be 
u = 0.25. The stress pulse is initially negative and then changes sign. The divergences at times 
corresponding to the longitudinal and Rayleigh wave arrival times have been cut off. 
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FIG. 3. Stress transfer kernels, J, for various models are shown as a function of time at a fixed 
distance z = 1. The thin line is the kernel for the quasistatic model and the squares that for the 
sound-travel-time delayed monotonic model. The open circles represent the scalar elastic model 
with To = 0.01 while the thick line represents the sharp pulse model with a = 0.5 and 7 = 1.5. For 
all except the quasistatic model, the stress transfer is zero for t < z. 
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FIG. 4. Series of avalanches in the quasistatic model in a system of size 256, as the load is 
gradually increased by just enough at each step to trigger the most weakly pinned site. The system 
is then allowed to evolve until motion stops before the load is increased again. The configuration 
of the crack front is shown by a thin line at the begining and by a thick line at the end of each 
avalanche to demarcate the sites which have moved. The position of the crack front in the figure is 
displaced vertically by a constant factor after each avalanche to differentiate between the individual 
avalanches. The initial almost straight configuration of the crack front at zero load is also shown. 
The avalanches shown occur at fractional loads (indicated on the right) in the range from 0.86 of 
the threshold load to the threshold load at which point the crack front starts moving. 



55 



o 
o 

'si- 




LJ 1 Q 

o o o o 

in o IT) 



FIG. 5. Space-time plot of a large avalanche in the quasistatic model in a system of size 512, 



with the points on the interface which moved ahead at each instant of time indicated. 
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FIG. 6. As the load is increased from zero towards the critical load, the avalanches 
that occur in the quasistatic model, are binned according to their size, which is defined 
as the number of distinct sites that move during an avalanche. The bin is defined by 
(bin size)/-\/2 < avalanche size < -\/2 (bin size). The number of avalanches in each bin are plotted 
versus the bin size for systems of size 1024. The slope of the log-log plot is 2.14zb0.3. 
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FIG. 7. Log-log plot of the mean duration of an avalanche versus the bin size is shown for the 
quasistatic model with system size 1024. The avalanches are binned according to their size which 
is determined by the number of sites that moved during that avalanche. The nth bin is defined by 
(bin size)/-\/2 < avalanche size < \/2 (bin size). The slope of the graph yields a dynamic critical 
exponent oi z = 0.74 zt 0.03. 
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FIG. 8. Log-log plot of the power spectrum versus the wave vector is shown for systems at the 
quasistatic critical load, ranging in size from 4 to 4096. They have been shifted along the vertical 
axis for clarity. 
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FIG. 9. Log-log plot of the power spectrum for system size 4096, averaged over 1000 samples 
and a linear fit over 0.004 <k < 0.35 which leads to 2C + 1 = 1.68 ± 0.04 
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FIG. 10. Log-log plot of the variance of the quasistatic threshold load as a function of the 
system size in the range from 4 to 8096. From a linear fit over the full range, f = 1.80 zt 0.05 while 
a fit for sizes from 256 to 8096 yields v = 1.72 ± 0.12 
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FIG. 12. Log-log plot of the mean velocity in the quasistatic model versus the the excess of the 
load above the critical value. The slope yields P = 0.68 ± 0.06. 
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FIG. 13. Space-time plot of an avalanche for (a) the time-delayed monotonic model (top) and 
(b) for the quasistatic model(bottom) starting from the same initial configuration of the crack 
front. The sets of jumps in the two figures are identical, but they occur at different times 
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FIG. 14. Velocity versus load for the time delayed monotonic kernel in a system of size 64, and 
a linear fit to the data. 
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FIG. 15. Space-time plot of an avalanche for the time-delayed monotonia model (open circles) 
and the stress pulse model with a = 0.5, 7 = 1.5 (with dots) from identical initial conditions. It 
can be seen that a large number of additional sites jump for the model with stress overshoots. 
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FIG. 17. Schematic of situation used in estimating the effects of a single small stress pulse. 
The space-time point A at which we estimate the probability of an additional jump due to a 
stress pulse from its neighbor jumping at a space-time point P, and P's "sound cone" showing P's 
isolation from the subavalanche, S, by a scale r, are all shown. 
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FIG. 18. The velocity of the crack front for the stress pulse model with a = 0.5, 7 = 1.5 is 
shown as a function of cycle number for a system of size 64 and extent 16 with periodic boundary 
conditions in the direction of motion. The crack front goes through the same random sample in 
each cycle and the velocity averaged over each cycle is measured while the external load is kept 
fixed. 
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FIG. 19. Log-log plot of the fractional decrease of the threshold force with scalar elastic stress 
transfer, (G^* — Gc(to))/G^*, versus tq is shown for systems of size 128. The slope of the linear fit 
is -1.56 
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FIG. 20. Mean velocity as a funtion of load for the scalar elastic model, in a single system of 
size 64, for various values of tq ranging from 1 to oo; tq = cxd corresponds to the monotonic model 
with time delays where we expect a continuous transition from the pinned to the moving phase. 
As can be seen, the finite systems exhibit hysteretic behaviour. 
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